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By adding minute concentrations of a high molecular weight polymer, liquid jets or 
bridges collapsing under the action of surface tension develop a characteristic shape of 
uniform threads connecting spherical fluid drops. In this paper, high-precision measure- 
ments of this beads-on-string structure are combined with a theoretical analysis of the 
limiting case of large polymer relaxation times, for which the evolution can be divided 
into two distinct regimes. This excludes the very late stages of the evolution, for which 
the polymers have become fully stretched. For times smaller than the polymer relaxation 
time, over which the beads-on-string structure develops, we give a simplified local de- 
scription, which still contains the full complexity of the problem. At times much larger 
than the relaxation time, we show that the solution consists of exponentially thinning 
threads connecting almost spherical drops. Both experiment and theoretical analysis of 
a one-dimensional model equation reveal a self-similar structure of the corner where a 
thread is attached to the neighbouring drops. 



1. Introduction 

Understanding the behaviour of polymeric free-surface flows is of enormous impor- 
tance for a wide variety of applications in the chemical processing, food and consumer 
products industries. Operations such as ink jet printing, spraying of fertilizers, paint- 
leveling, misting, bottle-filling and roll coating are all controlled by interactions between 
the non-Newtonian stresses in the bulk and capillary stresses at the deformable free- 
surface. Long-chained macromolecules are also ubiquitous in biological fluids, and very 
significantly affect the corresponding free-surface dynamics. If one places a small drop 
of saliva between two fingers and pulls them apart, the resulting liquid bridge does not 
collapse, but an extremely fine thread remains for several seconds. The lifetime of this 
bridge is intimately connected with the molecular weight and conformation of the pro- 
teins (mucins) and hormones in the saliva. Measurement of such filament lifetimes in 



2 



Ch. Clasen et al. 



biofluids such as mucus or saliva can be used as a fertility indicator Kopito & Kosasky 
(1979). 

A number of recent studies have promulgated the idea of using the capillarity-induced 

thinning of a liquid filament as a rhcomctric device for quantifying the properties of 
complex fluids in predominantly extensional flows (Bazilevsky et al. (1990), Stelter et 
al. (2000), Tripathi et al. (2000)). A typical configuration is shown in Figure 1. A liquid 
bridge of the fluid to be tested is initially formed between two coaxial cylindrical plates 
and a step uniaxial strain is then imposed on the bridge to extend it beyond the static 
(Plateau) stability limit. The liquid filament then undergoes a capillary-thinning process 
towards a final breakup event. The no-slip boundary condition at the endplates retards 
the radial flow near the end plates and thus imposes a well-defined initial axial perturba- 
tion or 'neck' on the liquid column which controls the location of the subsequent necking 
process. The time evolution in the local filament radius is monitored optically using ei- 
ther a laser micrometer or high-speed video imaging. To convert such measurements of 
filament evolution into an extensional viscosity or other material function characterizing 
free-surface extensional flows of complex fluids, it is necessary to understand the bal- 
ance of forces acting on the fluid filament. Capillary pressure drives the thinning process 
whilst viscous and elastic stresses (plus gravity) resist the necking of the fluid thread. The 
slcnderness of the fluid thread induced by the step strain means that a one-dimensional 
approximation to the equations can be useful and the large viscosity of many polymeric 
fluid systems means that inertial effects can often be neglected over the majority of the 
capillary-thinning process. Entov & Hinch (1997) provide a detailed discussion of the evo- 
lution of a cylindrical thread of viscoelastic fluid undergoing capillary-driven breakup. A 
central result of both this work and earlier studies is that there can be a lengthy inter- 
mediate regime in which inertial, viscous and gravity forces are all negligible and elastic 
and capillary forces balance each. In this regime the local extension rate in the neck is 
constant and the radius of the filament decreases exponentially in time. Measurement 
of this rate of thinning enables a direct determination of the characteristic relaxation 
time of the viscoelastic fluid. Such observations have been found to be in quantitative 
agreement with data obtained in extensional rheometers (Anna &: McKinley (2001)). 

This local elasto-capillary balance has also been observed in other experimental con- 
figurations including breakup of forced polymeric jets (Christanti & Walkca' (2002)) 
and gravity-driven drop formation in viscoelastic polymer solutions (Amaroucliene et 
al. (2001), Cooper-White et al. (2002)). In these experiments, the fluids typically have 
low viscosities and the initial stages of the process exhibit characteristics of a self-similar 
solution resulting from a balance between capillary pressure and fluid inertia (Day et 
al. (1998)) including a power-law decrease in the minimum radius. The rapid growth in 
viscoelastic stresses within the fluid neck leads to a cross-over to an exponential decrease 
in the minimum filament radius with time which is a hallmark of the elasto-capillary 
balance. However additional details of the self-similar solution (such as the spatial pro- 
file of the filament) are presently unknown. By contrast, in the absence of polymeric 
additives, the cross-over to a self-similar 'universal solution' balancing viscosity, inertia 
and capillarity are now well-explored (Eggers (1997), Rothert et al. (2001), Chen et al. 
(2002)). 

On close examination of a thinning viscoelastic jet, a string of tiny droplets can often 
be distinguished. This 'beads on a string' phenomenon was first described in Goldin et 
al. (1969), and has been reproduced in numerical simulations by Bousfield et al. (1986). 
A representative image of this viscoelastic jet break-up process is shown in 2. The jet 
consists of a series of cylindrical ligaments connecting spherical beads. As the jet is 
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Figure 1. Experimental images of a collapsing liquid bridge of polymer solution in a viscous sol- 
vent (plate radius Tio = 3mm, distance 13.8 mm). The surface tension is 7 = 37mN/m, the den- 
sity p = W26kg/m^. The solvent and polymeric contributions to the viscosity are ijs = 65.2Pas 
and rjp = 9.8Pas, respectively, the polymer timescale is A = 8.1s. Relative to the capillary 
timescale r = \/ pRnH this results in a Deborah number of De = 296. 
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Figure 2. High speed videoimage of a jet of dilute (0.01 wt%) aqueous polyacryamide solution 
(surface tension 7 — 62mN/ni) undergoing capillary thinning. The sharp-edged jet orifice is at 
the left of the image (radius Ro = 0.30 mm) and the free jet velocity is 30 cm/s. The polymeric 
contribution to the viscosity is qp — 0.0119Pas, and the polymer timescale is found to be 
A = 0.012s. This corresponds to a Deborah number of De = 18.2. 




convected from left to right, fluid is forced by capillarity from the thinning ligaments 
into the spherical droplets 

Most analytical studies of this structure have been performed using simplifying as- 
sumptions about the slenderness of the liquid jet, see Yarin (1993) for a review. How- 
ever, despite a considerable number of studies (see e.g. Goren & Gottlieb (1982), Entov & 
Yarin (1984), Bousficld ct al. (1986), Forest ot al. (1990), Shipman ct al. (1991), Larson 
(1992), Renardy (1994), Renardy (1995), Chang at al. (1999)) a full analytical description 
of the beads-on-string phenomenon is still open, even in the context of one- dimensional 
models. 

In this paper we seek a self-similar solution that encompasses the elasto-capillary bal- 
ance documented in experimental observations in liquid bridges, pinching drops and thin- 
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ning jets. We follow the spirit of some of the earlier work by employing two simplifying 
assumptions: 

First, we will consider the simplest canonical model for a dilute polymer solution, the 
so-called Oldroyd-B model (e.g. Bird et al. (1987)). This model is most easily derived 
from kinetic theory by treating a dilute solution of polymer chains as a suspension of 
non-interacting dumbbells, i.e. as two point masses or 'beads' at which the viscous drag 
from the solvent is imparted, connected by a Hookcan spring. The beads are assumed to 
be convected by the flow without slip, so in the extensional flow of a pinching thread, 
the beads are pulled apart, stretching the spring. In return, the dumbbell imparts a 
stress on the fluid, which results in a polymeric contribution to the stress tensor which 
is proportional to the force in the spring. The simplifying assumptions of a Hookean 
dumbbell are that the polymer chain exhibits a single relaxation time A, and that the 
spring is infinitely stretchable. The one-dimensional analysis of Entov & Hinch (1997) 
considered the more general case of a spectrum of relaxation times (corresponding to a 
non- interacting suspension of dumbbells with different spring constants). Their analysis 
showed that the rate of stretching in the liquid thread is governed by the dumbbells with 
the longest time constant and that all other modes relax and do not contribute to the 
dominant balance, so this approximation of a single relaxation time is not considered to 
be too limiting. Furthermore, a number of recent experimental studies (e.g. Spiegelberg 
et al. (1996), Amarouchene et al. (2001), Anna & McKinley (2001)) have utilized model 
dilute polymer solutions which are indeed very well described by a single time scale over 
a wide range of extensions. The additional assumption of infinite extensibility is bound 
to break down even for very long polymers as the trajectories of the two beads diverge 
exponentially in an extensional flow. In fact, it has been shown by Renardy (1994) for the 
model to be treated in this paper (and neglecting inertia) that a thread can never break 
up in finite time. In our study, we therefore disregard the final stages of breakup where 
the finite length of the polymers begins to affect the necking process. This regime has 
been considered by Entov & Hinch (1997), and the filament radius ultimately decreases 
linearly in time. 

The second simplifying assumption is that we are treating the flow inside the fluid 
thread as effectively one-dimensional (e.g. Forest et al. (1990)). This is consistent as long 
as the shape of the liquid column remains slender, i.e. the characteristic radial variations 
are small compared to the variation in the axial direction. This assumption is problematic 
near the ends of the fluid drops in the beads-on-string structure. Following Eggers & 
Dupont (1994), we hope to at least partially deal with this problem by keeping the full 
expression for the mean curvature in the Laplace pressure, which drives the breakup. 
This makes spherical drops exact static solutions of the equations, and ensures that at 
least the surface tension terms are correctly accounted for. 

We are left with a model that treats the liquid column as a set of one-dimensional 
continuum equations for the fluid flow coupled with equations describing the state of 
stress of the polymer chains in solution. A typical experimental situation would be that 
of a jet ejected from a nozzle, or a liquid bridge held between two circular end-plates. In 
all of the following, we will choose the initial bridge or jet radius i?o as a unit of length, 

1 /2 

and the corresponding capillary time t = (p_Rq/7) as unit of time, where 7 is the 
surface tension and p the density of the fluid. If Rq ~ 1mm, r is about 4 ms for a water- 
based solvent. Note that for a high viscosity fluid (also treated in this paper) other time 
scales such as the viscous scale r,; = rjR^y arising from a balance of surface tension and 
viscosity, might be more appropriate. However, to avoid confusion we will consistently 
use the inertial- capillary time scale. There still remain three independent dimensionless 
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parameters in the problem. The time scale A of the polymer is conventionally called 
a Deborah number, De, when made dimensionless with the characteristic time scale 
of the system. In the present study we thus have De = A/r. Note that the Deborah 
number is 'intrinsic' to the fluid thread because it is defined entirely in terms of material 
and geometric parameters. It does not contain the rate of stretching in the fluid, since 
the flow is not forced but is frc;e to select its own rate of stretching, which may be 
spatially and/or temporally inhomogeneous. The other two dimensionless parameters 
represent the relative contributions of viscous stresses from the solvent and the polymer. 
There are a number of possible representations for these parameters. The total dynamical 
viscosity for a dilute polymer solution characterized by the Oldroyd-B model is given by 
T]o = Vs + Vp and the relative importance of viscous effects can thus be characterized by 
the Ohnesorge number Oh = qo /^p'fRo and the solvent viscosity ratio S = r?s/r?o- An 
alternative representation which we use below is to separate the relative dimensionless 
contributions of the kinematic viscosity Vs = OhS of the solvent, and the polymeric 
contribution to the viscosity Up = Oh{l — S). All these material constants have been 
made dimensionless using Rq, 7 and r as the characteristic scales. 

Even if one only considers the case of a fixed initial condition, and neglects gravity 
(say in a liquid bridge configuration), the problem is still one of daunting complexity. In 
particular, the behaviour of the liquid bridge need not be spatially uniform throughout: 
in some regions the rate v' at which fiuid elements are stretched may be smaller than 
l/{2De). In this case polymer relaxation overcomes flow stretching, and the polymers 
remain coiled, leading to Newtonian-like behaviour. In other parts of the flow stretching 
overpowers relaxation, and the polymers become stretched in unison with the exponen- 
tially diverging paths of fluid elements. In these regions, the polymeric contribution to 
the stress grows accordingly, and threads form. The transition from a flow regime with 
2v'De < 1 to one with 2v'De > 1, where exponential stretching of polymers may occur, 
has been called the "coil-stretch" transition by dcGennes (1974). 

It is very difficult to describe this complexity in full generality, so in the following we 
will restrict ourselves to the case of large De » 1, implying that the non-Newtonian 
polymer contribution is significant at all times. Physically this means that De is much 
larger than the initial time scale of the liquid bridge's evolution, which is set by the linear 
stability of the fluid thread. At low viscosities. Oh = Pg + < ^, this time scale is 0(1) 
by virtue of the chosen time scale for non-dimensionalization. By contrast, for fluids with 
large viscosities it is set by Tj, = r/Oh, and we thus require De/Oh ^ 1. Note that with 
the present scaling, 1/Oh^ is the ratio of the external length scale Rq and the intrinsic 
scale of the fluid, £^ = iy^/{-fp). 

With these assumptions, we can divide the transient evolution of the jet into distinct 
stages as shown schematically in Figure 3. The dimensional time i elapsed is scaled with 
both the polymer relaxation time A and the appropriate time scale for evolution of the 
fluid thread Tch (which depends on the value of the Ohnesorge number). Fluid elements 
in the jet evolve along pathlines of constant slope X/tch = De^h- In each regime, we 
can provide a simplified description by balancing the dominant terms in the governing 
equations. The focus of this paper will be on the sector characterized by large Dcch- 

In the early elastic time regime, t <C 1 ^ De, there is no significant decay of polymer 
stretching. The fluid thus responds as a neo-Hookcan elastic material. This allows the 
effect of the polymers to be written as a local contribution to the pressure, given in 
terms of the interface shape (i.e. the local accumulated strain) alone. The parameter 
determining the magnitude of this contribution is the dimensionless elastic modulus 
of the material G = Vp/De, which is (up to universal constants) proportional to the 
polymer concentration. Depending on the viscosity, the dynamics of the bridge can be 
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quite complex. In particular, for low viscosities, capillary waves can travel along threads 
and rebound off drops (Li & Fontelos (2003)). Threads are also shown to support elastic 
waves. 

For De » t » 1, as polymers become sufficiently stretched to counter surface tension 

forces, the simplified, local system of equations converges to a stationary solution, main- 
tained by the stress in the polymers with no possibility of relaxation. This stationary 
solution, originally found by Entov & Yarin (1984), already exhibits the beads-on-string 
structure, but with a thread of radius hi^read = (0/2)^/^ to be computed in section 3.2. 
The transition of the initial evolution to the region marked "quasi-static" in Fig. 3 thus 
occurs approximately when this radius is reached. The name "quasi-static" refers to the 
fact that the solution can only be regarded as stationary on time scales much smaller 
than the polymer relaxation time. 

Indeed, to proceed beyond this stage one has to take the viscoelastic relaxation of the 
polymer chains into account. The structure of the solution is that of cylindrical filaments 
which thin at an exponential rate cxp(— t/3-De) as a result of the local balance between 
elasticity and capillarity. The filaments connect an arbitrary distribution of droplets, 
which approach a static, spherical shape. A similarity solution describes the crossover 
between the cylindrical thread and the neighbouring droplet. Towards the thread, the 
solution asymptotes towards a constant thickness, in the direction of the drop it fits 
onto the spherical shape of the drop. Since the asymptotic shape is that of threads of 
vanishing radius connecting perfect circles, a corner develops, with the interface slope 
diverging on the side of the drop. Thus the description of this part of the flow is beyond 
the capabilities of a slender-jet model. However, by keeping the full curvature term, we 
still hope to obtain a reasonable description of this transition region. At low values of 
De, the initial stages of the necking process at short times t < De < 1 are controlled by 
the linear stability of a viscoelastic fluid thread (Middleman (1965), Punada & Joseph 
(2003)); however ultimately at longer times De < 1 < t there will still be a crossover to 
the exponential necking regime, provided the polymer chains are sufficiently extensible. 
We do not consider these regimes further in the present study. 

Our paper is organised as follows: In the next section we develop and motivate the 
lubrication equations to be used for the remainder of this paper. A numerical simulation 
illustrates the regimes to be analysed below. The third section is devoted to the study 
of the neo-Hookean regime, where polymer relaxation can be neglected. First we derive 
a simplified, local description from the model equations. This local description is then 
used to compute the asymptotic thread radius, and to describe the propagation of elastic 
waves on the thread. The fourth section deals with the long time regime at finite De for 
which exponential thinning of threads is observed. After giving a qualitative description 
of the shape and fiow inside the thread, we introduce a similarity description valid in the 
corner where a cylindrical thread meets a spherical drop. If De is large enough to make 
only clastic and surface contributions relevant, we can compute all but one; of the free 
parameters of the solution. This last parameter, the thread radius, can be estimated by 
matching to the early-time regime. The numerical results are compared with experimental 
observations using a dilute solution of monodisperse polystyrene which is well-described 
by the Oldroyd-B constitutive model. Measurements of the evolution in the mid-filament 
radius and the evolution of the spatial profile of the filament are well-described by the 
theory. In the final section we discuss work that remains to be done within the framework 
of the present model, as well as perspectives for inclusion of other effects that lie beyond 
it. 
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Figure 3. Schematic diagram showing the transient evolution of viscoelastic fluid threads. 
Trajectories correspond to rays of constant slope, given by the the 'intrinsic' Deborah number 
De — X/tch, where Tch is the characteristic capillary timescale (depending on whether viscous 
or inertial effects are more important in the bulk). The focus of this paper is on the sector 
corresponding to large slopes, and more details are given in the text. For small slopes the 
behaviour is Newtonian for a significant portion of the evolution, until the stretch rate v' exceeds 
the inverse polymer relaxation time. 



2. Model and simulation 

2.1. One- dimensional equations 

In this paper we confine ourselves to the study of a simphfied version of the Oldroyd B 
model for polymeric liquids, assuming that the radius h{z,t) of the liquid column varies 
slowly. Thus the variation of hydrodynamic variables inside the column is also small, and 
we can confine ourselves to the leading order approximation in an expansion in the radius. 
For example, v{z,t) below is the axial velocity at the centre of the jet. A derivation of 
the relevant equations has been given in Forest et al. (1990), so we just give the final 
result and briefly discuss its physical significance: 



f + |W=0 (2.1) 

dv dv Ok ^ 1 d f,^dv\ 1 d 
dt dz dz ^ h-^ dz 

^ + " 97 ^ ^9^"- ^^D-ed-z-We ^''^^ 
Equation (2.1) expresses volume conservation, (2.2) is the momentum balance equation 
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in the one-dimensional approximation. The first term on the right of (2.2) is the gradient 
of the Laplace pressure, which is the main driving force. With the present scaling, the 
capillary pressure gradient is simply the gradient of the mean curvature k, which for later 
convenience will also be written as 



The second term on the right of (2.2) is the Newtonian contribution to the viscosity, 
multiplied by the dimensionless viscosity Vg of the solvent. Finally the last term is the 
polymeric contribution; azz and the diagonal terms of the extra stress tensor. 

The other components of the polymeric stress tensor do not enter at leading order. Using 
(2.5), equation (2.2) can finally be rewritten such that the inertial terms on the left are 
balanced by gradients of the tension in the thread: 



For very viscous fluids the inertial terms on the left are negligible, and the terms 
in parentheses equal a (generally time-dependent) constant. If velocity gradients dv/dz 
are small, (2.3) and (2.4) describe an additional Newtonian contribution Vp to the total 
steady-state shear viscosity u = Vg + ^p- The presence of polymers also results in fluid 
viscoclasticity and the polymeric stress in the fluid responds with a dimensionless time 
delay of De, which is represented by the relaxation terms UijjDe on the right of (2.3) 
and (2.4). Finally a crucial term for the physics of the following is the first term on the 
right of (2.3) and (2.4), which describes the interaction of the polymer with the flow. 
In an extensional flow dv/dz is positive, so the stress in the axial direction grows as 
the dumbbells modeling the polymeric contribution to the stress are stretched, while it 
decays in the radial direction. 



We are now in a position to study the behaviour of the model for various initial conditions 
and to compare to experiment. First, we simulate the decay of a long, initially unstretched 
cylinder of fluid. For our simulations we have used a numerical code analogous to the 
one developed earlier by Eggcrs & Dupont (1994), Eggers (1997) for Newtonian flows. It 
is fully implicit and uses adaptive regridding of the numerical mesh to fully resolve the 
fine structure of the fiow. This is crucial to be able to describe some of the last stages of 
thread formation to be investigated in detail in section 4. We found that the demands 
on the solution of the implicit equations are much greater than in the Newtonian case, 
owing to a larger range of time scales in the fiow. In general, several iterations of a 
Newton scheme were necessary for convergence, and significant restrictions had to be 
put on the time step. To further test for possible problems inherent in our numerical 
scheme. Another explicit code was also developed independently. In addition, all fields 
were represented in a uniform grid, as opposed to the staggered grid of the implicit code. 
The explicit code performed quite well except for the highest viscosities, where the time 
step imposed by the Courant condition became prohibitively small. In all our tests, we 
found no significant disagreement between the results of the two codes, so below we will 
no longer specify which one was used to obtain a specific result. 

Figures 4 and 5 give an idea of the typical behaviour of a liquid filament in the absence 
of gravity, described by (2.1-2.4). The parameters were chosen to be identical to those of 
Fig. 2 in Chang at al. (1999); that is in the present scaling, Vg = 0.79, Vp = 2.37, and De = 
94.9, except that our initial perturbation is smaller. The periodic boundary condition is 




(2.5) 




(2.6) 



2.2. Beads on a string 
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Figure 4. A typical series of profiles with periodic boundary conditions and period L = 4-w. 
The dimensionless parameters are Vs = 0.79, Vp = 2.37, and De = 94.9. A very small sinusoidal 
perturbation was added in order to make the filament collapse. After the rapid formation of 

the beads-on-string structure, one observes the slow thinning of the thread. The relative di- 
mensionless times of each profile are (a) 0.0, (b) 31.6, (c) 158.1, (d) 316.2, (e) 632.5 and (f) 
948.7. 
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Figure 5. The minimum radius corresponding to the profiles shown in Fig. 4. One clearly ob- 
serves a rapid initial motion for times small compared to De = 94.9, followed by an exponential 
thinning at rate l/{3De). 

used here with the period L = An. The relaxation time scale for the polymer solute (A) is 
approximately 95 times the capillary time scale, so at short times, viscoelastic relaxation 
effects can be neglected. Additionally, the Ohnesorge number is Oh = Vs + Up = 3.16, 
so the flow is significantly damped by viscosity. Fig. 4 shows a sequence of profiles for a 
length of two periods (Stt) at different times. After an initial period described by linear 
theory, a thin thread is formed between two round drops. The filament necks in the 
central rc;gion and th(^ fluid expelled from this neck is accumulated in the two beads. The 
evolution in the minimum radius of the thread with time is shown in Fig. 5. There is 
an initial sharp decrease, until a sudden change in the slope occurs. This happens when 
the elastic stresses which build up in the deforming liquid bridge become dominant and 
a local elasto-capillary balance leads to the formation of a thin filament. The rest of the 
evolution consists of an exponential thinning of the filament, best seen in the logarithmic 
plot of Fig. 5. The slope drawn into the figure corresponds to the theoretical prediction 
of section 4. 

Structures very similar to Figure 4 are shown in the experiment of Figure 2 and 
have also been observed by Bazilevskii ct al. (1981), Christanti & Walker (2001) for the 
decay of a liquid jet of polymer solution ejected from a nozzle. Good agreement with 
numerical simulations of one-dimensional models very similar to ours have been reported 
by Yarin (1993). However, it is very difficult in practice to produce liquid cylinders 
without stretching the polymers, since there is considerable shear inside the capillary tube 
and the nozzle. On the other hand, the shear flow inside the capillary is very difficult to 
model in the framework of the present one- dimensional description. Furthermore, each 
bead and ligament shown in Fig. 2 corresponds to a 'snapshot' at a different elapsed 
convective time At = Lperiod/vjet Therefore, to compare quantitatively to experiments. 
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Figure 6. A comparison between experimental profiles (left), obtained from digitising the im- 
ages of Fig. 1 and corresponding simulations (right), including the initial stretching and gravity 
(acting downward). The time difference between two consecutive images is 10 sec. The re- 
duced variables corresponding to the experimental parameters are Vs = 193.2, z/p = 29.04, and 
De = 296. Note that there are more profiles from simulation to reach approximately the same 
thread radius, corresponding to a slight overestimation of the experimental time scale. 



we prefer to use a setup that allows for a more quantitative description of the stretching 
history of the fluid column. 

2.3. Liquid bridge 

We now turn to the capillary thinning and breakup of a liquid bridge that has been sub- 
ject to a very rapid initial stretching. During this extensional step strain process, which 
lasts 0.05 sec, hardly any polymer relaxation takes place, so the exact manner in which 
the the plates are pulled apart is not very important. In fact, the initial stretching is well 
described by a simple model neglecting any spatial structure (Anna & McKinley (2001)). 
However, the simulation within the one-dimensional model is somewhat subtle, owing to 
a difficulty in imposing the boundary conditions at the end-plate. Namely, the no-slip 
boundary condition enforces a vanishing tangential velocity at the end-plate, and conse- 
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Figure 7. The logarithm of the normahsod minimum radius corresponding to the experimental 
and theoretical profiles of Fig.6. The thick line is the theoretical prediction for the slope — 1/(3A). 

quently dv/dz = 0, while this is not true in the one-dimensional model. Since the stretch 
rate dv/d/. is in fact large over most of the bridge, this creates a thin boundary layer of 
the full three-dimensional flow near the end-plates. Failure to correctly implement this 
boundary layer leads to a detachment of the interface from the ends within the lubrica- 
tion model. Following Stokes ct al. (2000) , we have avoided this problem by introducing a 
supplementary viscosity which strongly increases near the ends. This position-dependent 
viscosity has been constructed by matching to a three-dimensional squeeze flow near a 
solid wall. The effective "freezing" of the fluid prevents any lateral slip along the bounding 
wall. 

Figures 6 and 7 allow for a direct comparison between simulation and experiment. 
The first digitised profile is taken just after cessation of stretching, after that a profile 
is shown every 10 sec. Theory and experiment show good agreement in all the basic 
features of the flow, such as the sagging under gravity and the formation of the thread. 
Two subtle differences can be seen: first, the timcscalc of the simulation is off by about 
20%, so 9 experimental profiles are shown but only 7 theoretical ones, at which point 
about the same minimum thread radius is reached. This discrepancy, also seen in Figure 
7, is quite acceptable considering that no adjustable parameters were introduced. Most 
of the difference stems from the early-time development, while the asymptotic slope of 
both simulation and experiment agrees well with the theoretical prediction. 

The behaviour of hmin at early times is markedly different from that of a free jet, as 
discussed in detail in Anna & McKinley (2001). The rapid early decrease of hmin seen 
in Figure 5 is absent, since stresses are already large owing to the initial stretching of 
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the liquid bridge. On the contrary, some of this initial stress has to relax before further 
thinning can start, as seen in the plateau for the earliest times. Evidently there are 
some subtle features of the experimental relaxation processes which are not modelled 
correctly by our single-mode Oldroyd-B model. The initial stretch is also responsible 
for the absence of drops ("beads") in the middle of the thread that formed on the free 
jet, cf. Figure 4. The reason is that the initial stretch is uniform, and this uniformity is 
conserved by the exponential stretching regime. 

The second difference between the experimental and theoretical profiles of Figure 6 is 
that at the same minimum thread radius the corner between the thread and drops at the 
end is sharper. We will return to this when we discuss the structure of the corner region 
in detail in section 4.2. 



3. Early time asymptotics 

3.1. Local description 

In this section we show that for early times < <C 1 <C De the non-Newtonian contribution 
to the stress behaves like that of a neo-Hookean elastic solid, and can be expressed 
directly through the profile shape h{z). This proves the phenomcnological argument of 
Entov & Yarin (1984), in addition to identifying the correct value of the elastic modulus 
of the material: 

G = Vp/De. (3.1) 

With the present scaling this is the same as the elastocapillary number introduced in 
Anna & McKinley (2001) and can also be written G = Oh{l — S)/De. Following the 
procedure in Fontelos (2003) we rewrite (2. 3), (2. 4) in Lagrangian coordinates. This idea, 
introduced in Renardy (1994), is useful since the stretch of a fluid element can be written 
as Za = (Renardy (1994), Eggers (1997)). This Lagrangian stretch u = Za is the 
derivative of the particle position z with respect to the label a of the particle in a reference 
configuration. Using the Lagrangian definition of the fluid velocity Zt = v one finds that 
Vz = Ut/u and thus 

a^z + De -—— - 2—azz = 2vp— , (3.2) 
\ at u ) u 

( darr Ut \ Ut , , 

CTrr + Ue I —-— H CTrr = -^p— , (3.3) 



which can be rewritten as 



dt u 



<_9/l_!_ \ Ut 

' '"dt[-^''"^-r^''u^ 



dt 



& '^-TT. [ue^UrT^ = —Gut ■ 



If we assume that the initial stresses are zero (otherwise they remain as initial condi- 
tions) we can integrate the equations to find 

azz = '^Gv?e~^ I e^^dr = -Gu^e'^ [ e"^ f ^ | dr , 
Jo Jo WJr 

Jo 



arr = —Gu I e'^Urdr 
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Integration by parts yields 




G 2 —L- f 

ut De Jo 



—dr , 



u- 



(3.4) 



C7rr = — Lr + G t-C Oe -|- -—U 6 




f e udr , 





(3.5) 



where the integrals are taken along Lagrangian paths, i.e. for constant particle label a. 
However, for t ^ De the integral can be neglected relative to De and the exponential 
is close to unity. For the case that the initial deformation of the bridge is small, i.e. 
that Uq is almost constant, (3.4,3.5) only contains Lagrangian labels at the same time t, 
which can simply be replaced by ordinary Eulerian coordinates by applying the inverse 
transformation. Remembering also that with the scalings used throughout this paper 
uo = Ro = 1, we find 



This expression is the well-known result for the neo-Hookean rubber material if we rec- 
ognize that the Lagrangian stretch isu = l/h"^. Thus finally the equation for the velocity 
is 



where K is defined in (2.5). Of course, (3.7) can also be confirmed directly by applying 
the limit t <^ De to the Eulerian equations (2.3,2.4), but it is hard to guess a priori 
without the use of Lagrangian coordinates. 

As a numerical test of the quality of the local approximation, we performed two sim- 
ulations similar to that of Fig. 4, but for two different Deborah numbers. Fig. 8 shows 
the evolution of the interface profile at different times for De = 94.9, Up = 2.37 (left 
column), and De = 9490, i/p = 237 (right column). The full lines are the solution of 
equations (2.1-2.4), while the dashed ones were obtained by replacing (2.2) by (3.7). For 
the simulation corresponding to moderate Deborah number, the duration t is of same 
order of the relaxation time De. And as we expected, the solution of the full equations 
and the one of the local approximation differ for times t > De. On the other hand, for the 
case of large Deborah number, t <C -De, the agreement of the two solutions is excellent. 
Indeed, we can hardly distinguish the dashed lines from the full ones in the right column 
with our naked eyes. 

Throughout the process of filament formation, we find almost perfect agreement with 
the local model provided De ;» 1. Thus the local approximation is an extremely useful 
tool to investigate the early-time dynamics t < 1 < De and the formation of the basic 
bcads-on-string structure. Wc will therefore study the stationary filament solutions of the 
local model in the following subsection. After the filament has formed, the neo-Hookean 
elastic response of the local model leads to a stationary profile, while the efii'ects of fluid 
viscoelasticity captured by the full equations leads to the filament continuing to thin 
exponentially in time. 

At low viscosities, Oh < 1, the initial evolution can be considerably more complicated, 
but is still fully described by the local equations if De is sufficic^ntly large. This is due to 
inertial effects also being important, so as the thread is formed, fluid may rebound from 
the drops, and an additional 'secondary' drop forms in the middle of the filament. This is 
also seen in experimental observations of jet breakup (Christanti & Walker (2001)) and 
drop pinch-off (Cooper- White et al. (2002)). In addition, capillary waves may develop, 
which we describe in subsection 3.3. 



h^ {(Tzz - (Trr) = G (l/Zi^ - h^) + 0{t/De). 



(3.6) 




(3.7) 
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Figure 8. Comparison of interface profiles with period L = 20 and Vs = 0.79 between the 
full equations (full) and the local approximation (dashed lines): left column De = 94.9 and 
Vp = 2.37, and right column De = 9490 and Up = 237. The relative times for each profile are 
(a) 31.6, (b) 79.1, (c) 158.1, (d) 316.2 and (e) 474.3. Evidently, for t < De the agreement is 
excellent (right column). 



3.2. Static solutions 

In the neo-Hookcan limit considered here, De oo, the elastic stresses never relax, so 
at long times t 3> De » 1, surface tension is balanced by permanent elastic stresses to 
form a stationary solution, which obeys the equation 

h^K + G{l/h^ - h'^) = C, (3.8) 

which comes from integrating (3.7) once while dropping inertial terms. Apart from the 
appropriate boundary conditions on h, (3.8) has to be solved with the constraint of 
volume conservation 

TT [ h^ = V, (3.9) 
which is needed to find all the constants of integration as well as C. 
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In particular, (3.8) allows for solutions with constant h 

h + G{l/h'^-h'^) =C, (3.10) 

where we have assumed a constant initial height Hq = 1. The solutions of equation 
(3.10) correspond to the thin cylindrical threads of constant radius described earlier. To 
determine their radius, one needs the value of the constant C, which only follows from 
a solution of the full system (3.8) with proper constraint (3.9). Entov & Yarin (1984) 
have described a procedure which matches solutions of (3.10) to almost circular drops in 
the limit of small elastic modulus G, from which the thread radius can be determined 
analytically. Below we will reproduce their calculation, which will turn out be relevant 
for the calculation of thread thinning as well. 

First, in the drop region the radius h is of order one and one can neglect the terms 
multiplied by G in the limit we are interested in. Thus one is left with the contribution 
from surface tension alone, and the solution is a spherical drop. 



h{z) = R^l - (1 + z/i?)2. (3.11) 

Here R is the radius of the drop, which is set by the volume constraint (3.9), and must 
be a quantity of order one. The other constant of integration was used to place the edge 
of the drop, which merges onto the thread, at the origin. In the following we will assume 
that the filament is to the right of the drop, in the region of positive z, while the matching 
takes place in the neighbourhood of the origin. 

To match the drop (3.11) onto a constant solution, we have to analyse (3.8) in the limit 
G 0, guessing that in the crossover region h is still small. Thus h'^ can be neglected 
relative to l/h"^ and we can eliminate the elastic modulus from the equation with the 
scaling 

h = G^/^H, z = G^/^Z, C = G^/3Ci, (3.12) 
where now all rescaled quantities arc of order unity and H obeys the equation 

H^K{H} + 1/H'^ = Ci. (3.13) 

Denoting derivatives with respect to X by primes and multiplying (3.13) by H' /H^ one 
finds the following first integral. 



H' = 



{Ci+2m/Ri - l/(2ff2))2 
where Ri is a constant of integration. Here we have used 



1/2 

(3.14) 



1 \ Hz H^Hzz 



The the constant Ci and Ri in the equation (3.14) have to be chosen so that the solution 
matches both the drop and the thread. 

We begin with the matching to the drop for ^; < 0. For large Ri there is region with 

rY^ ^H^Ri (3.15) 

in which (3.14) is approximated by 

H' w -Ri/H 



with solution H = y/2Ri\z\. By comparing with the expansion of the drop profile (3.11) 
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for small |^| we find 

Ri = R/G^/^, (3.16) 

so Ri is indeed large for small P. In terms of the variable z, the matching region (3.15) 
is defined by 

1< |z| < (3.17) 
On the other (thread) side of the problem (z > 0) the slope has to go to zero and thus 

C, + 2Hl/Rl - l/{2Hl) = 2i7o, 

giving Ci in terms of Hq. To finally find _ffo: one has to make sure that small perturbations 
around the constant solution do not grow. Putting H — H,^{1 + 5) , we find 



5' - [(1 - 2Hq/Rx - l/{2Hl))5 + 0(<52)] 



1/2 



For 5 not to grow the expression in parentheses must vanish, and since Ri becomes large 
in the Umit G ^ 0, we are finally left with Hq = 1/2^3 and Ci = 3/2 V3. Thus the 
central result of our analysis is that the thickness of the thread becomes 

f G\ 

hthread = [-:^ \ ■ (3.18) 

For the parameters of Fig. 8 we find h-ihread = 0.232 from (3.18), to be compared with the 
observed value of hmin = 0.253. The quality of the approximation quickly improves as 
the scale separation between thread thickness and drop size becomes even more complete. 

3.3. Linear stability 

The analyses of the previous two subsections can be used to perform a linear stability 
analysis of the stationary threads that form in the limit of very large De. This will be 
a good approximation on time scales much smaller than < < 1 <C De, over which the 
thread maintains an almost uniform axial profile. This is a relatively simple calculation 
to do, since the base solution is constant in both space and time. Thus we put 

h = ho-\-h, v = v, (3.19) 

where Hq is constant and h, v are small perturbations around the base solution. Then 
the linearised equations turn out to be: 

dh 1 , 9^7 

which can be combined to the following equation for h alone: 

Using a wave ansatz of the form h{z, t) = e expi{kz — cot) we find the dispersion relation 

'^^ = - + + + 2.iVshowe. (3.21) 

2/10 "o 

The first term on the right comes from surface tension and describes the usual Rayleigh- 
Plateau instability of a Newtonian liquid jet. If we disregard the polymer contribution 



18 



Ch. Clasen et al. 



we find complex u) for kho < 1 and thus exponential growth of perturbations. The 
presence of elastic stresses shifts the critical value of k to slightly smaller values of k, 
thus stabilising the flow. If however this critical value as well as the maximum growth rate 
are compared to a Newtonian fluid with the same zero shear rate viscosity 770, one finds 
that the instability grows faster (e.g. Middleman (1965), Goldin et al. (1969), Chang at 
al. (1999)). 

In the thin thread region described by (3.18), on the other hand, one is effectively 
dealing with an elastic medium and can see the propagation of elastic waves. Only taking 
into account the elastic contributions one finds 

-,1/2 



± 



14 



3iUshok/2. 



(3.22) 



Hence for small viscosities there is little damping and waves are very easily excited on 
the threads, which causes considerable numerical problems to properly resolve them. In 
the limit of zero viscosity the elastic wave speed is given by 



Vel = 



hi ■ 



(3.23) 



4. Late time asymptotics 

4.1. Thread thinning 

The formation of threads described in the previous section is a result of the interplay 
of surface tension and clastic forces. Any additional thinning of the thread results in an 
extensional flow, resisted by further stretching of the polymers. On times longer than 
De however the stress gradually relaxes, and the thread thins at an exponential rate 
p. This rate is easily determined from a balance of surface tension and elastic forces 
in (2.2), assuming a spatially constant profile (e.g. Bazilevskii et al. (1981), Renardy 
(1995), Entov & Hinch (1997)) 



h{z,t) = hoCxp{-(3t). (4.1) 

From volume conservation (2.1) one finds that the extension rate d^v = 2/3 in the 
thread is constant. The exponential growth of the axial stress is described by (2.3), 
and assuming a spatially constant one immediately finds 



azziz, t) = ao exp[(4/? - De)t]. (4.2) 
The radial stress arr decreases exponentially and does not figure in the balance. Re- 
membering that capillary pressure is balanced with azz in (2-2) and rewriting the pressure 
gradient according to (2.5) one finds 

h + h'^azz = C, (4.3) 
performing one spatial integration. For the balance (4.3) to be consistent azz must 
grow like l/h, and thus /3 = l/(3£'e), implying that the constant of integration C itself 
decays like 

C = aiexp[-i/3£>e]. (4.4) 
This means that the thinning rate of the thread is directly related to the time scale of 
the polymer, providing for a convenient experimental probe. 

To obtain a clearer physical picture, we plot in Figure 9 the thread radius, normalised 
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Figure 9. Profiles of the radius, the velocity, and the axial stress, normalised by their val- 
ues in the thread, in the regime of exponential thinning. The parameters are De = 3164 and 
I's = I'p = 95.6. The extension rate inside the thread is 2/{3De). 



by the exponential factor £ = exp[— t/3-De], the velocity in the thread, and iazz for 
a number of different times during the thinning. As predicted, h{z) and a{z) are very 
nearly constant over the thread, and they collapse nicely as anticipated by the above 
scaling laws. Furthermore, the extension rate dzV has the constant positive value 2/3 
inside the thread, expressing the fact that fluid is expelled from it. In response, the stress 
a^z grows to large positive values, so both contributions on the left-hand side of (4.3) are 
always positive. This means that the constant of integration in (4.3) (which below we will 
compute explicitly in the limit of large De) has to be kept for a consistent balance and is 
needed to determine the stress and thus the extensional viscosity from a measurement of 
the thread radius. For large De this makes the stress twice as large (apart from a difference 
in sign) as determined from the customary direct balance between capillary and elastic 
stresses (e.g. Bazilevskii et al. (1981), Entov & Hinch (1997), Anna & McKinley (2001)). 
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Since the total length L of the thread is a constant quantity, the maximum value of 
the velocity 

L , , 

Vmax ~ (4.5 

■iUe 

behaves like the inverse of the Deborah number, and is thus small in the limit that is the 
chief focus of this paper. A number of experiments (e.g. Bazilevsky et al. (1997), Anna & 
McKinley (2001)) have confirmed the prediction (4.1). In particular in Anna & McKinley 
(2001) the relaxation time was determined independently, and the thinning rate was 
found to conform with the prediction /3 = l/{3De). At both ends of the thread, the 
velocity and the stress fall to zero very sharply, and the profile merges onto a static 
drop with radius R. Next we will focus on this transition region, whose scale is set by 
i = exp{-t/3De). 

4.2. The corner region 
According to the scalings found in the previous subsection, which implies the existence 
of a small length scale i = exp(— t/3£)e), it is natural to look for solutions of (2.1-2.2) of 
the form 



h{z,t) = eh{z,t) 

v{z,t)=v{z,t) (4.6) 
a,,{z,t)=£-^7f,,{z.t). 

where z = i~^{z — zq). The origin zq must asymptotically lie in the similarity region. 
A convenient choice is the position of the extremum of the velocity in the limit £ — > 0. 
Since arr is exponentially small inside the thread, it can be left out of our analysis. Thus 
the equations for h,v, and azz are 



2ihh 



2£ 



5Dt 



I 



3De 



(^h — zhh ^ 



{azz + za',,)+v^ 



vh 



K{h} + {3i^sv' + azz) 



— (2fp'(;' - o-^z). 



(4.7) 



where the prime refers to differentiation with respect to the similarity variable z and 
the overdot indicates differentiation with respect to time. Towards the thread our scaling 
ensures that h,v, and a^z tend toward the constants ho,vo, and ctq, respectively, as 
J ^ oo. (Without loss of generality, we assume that the thread is to the right of the 
transition region.) The length scale £ becomes exponentially small for t oo and for any 
finite z, terms proportional to £ or £'^ may be dropped. Looking for time-independent 
solutions of (4.7) and integrating once, we find 



vh = vohg 

liK{h} + }i{3vsv' -\-azz) = hQ + KqCJo 



(4.8) 



Eliminating v and Uzz, we end up with 



h^K{h} + h^{-a2h'/h^ 



az/h ) = ai, 



(4.9) 




Figure 10. A schematic of the different regions near the corner. 



where ai = + hQ<To,a2 = Gi^gUo^oj ^-^id = ctq^o- general, equation (4.9) is no 
longer valid inside the drop, because here the profile h{z,t) is almost stationary, thus 

h varies like £^^, implying that the time derivative h is proportional to l/{£De). In 
addition, the similarity variable J is also of order l/£ inside the drop. Therefore we have 
to distinguish three different regions, as shown in schematically in Figure 10. 

(a) The similarity region near the thread described by (4.9). 

(c) The solution inside the almost spherical drop, where the velocity and stress fields 
relax towards zero. 

(b) A transition region connecting the dynamical region (a) with the static region (c). 
However in the limit of large De which is the main focus of this paper, velocities are of 

the order of l/De, and the contribution —02/1 /h in (4.9), which comes from the viscous 

stress, can be neglected. The same is true for the left hand side of the second equation 

(4.7), as we will show in more detail in the appendix. Physically this limit corresponds to 

a situation where the fluid flow becomes irrelevant to the problem, which is dominated 

by surface tension and elastic contributions. Since the solution of (4.9) matches onto the 

drop, h must become as large as desired for a sufficiently large z. But this means that 
2 2 

the contribution from h 'a^z = as/h is negligible compared to the constant ai outside of 
the similarity region. Thus although azz can no longer be written as a simple function of 
h in regions (b) and (c), it does not contribute to leading order outside of the similarity 
region (a). In summary, this implies that (4.9) with 02 = is valid everywhere in regions 
(a)-(c). 

The constant 03 can be scaled out by putting 

h{z, t) = al^^giC = z/al^^ t), (4.10) 

and so we finally find 

g'K{g} + 1// = Ci (4.11) 

where Ci = ai/a^^^. 

Equation (4.11) is identical to (3.13), describing the stationary thread in the local 
approximation. The matching described in section 3.2 

can thus be applied with very slight modifications. Multiplying (4.11) by the integrating 
factor g'/ g^ and integrating once, we obtain 
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g' 



' V 

.(Ci+2<7Vi?i(f)-l/(252))2 
We seek values of Ci and Ri{t) such that 



1/2 

(4.12) 



g{C,t) Const. whenC>l, (4.13) 
g{C,t)^{2R)i{-£-^a^^^\)i when C « -1, (4.14) 

implying matching with the filament and the drop respectively. In order to match to 

1 /3 

the stationary drop of radius R we put Ri{t) = R/{ia^' ) in analogy to (3.16), and in 
the limit of f — > oo we arrive at the equation 



9' = 



4g2 



The procedure to match with the filament described in section 3.2 establishes that 
Ci = 3/2V3 and 

(To = 2//io, ai = 3ho. (4.16) 
The transition region (b) between the drop and the thread is by virtue of (3.17) defined 

by 



1< 1^1 < R/e. (4.17) 



In this region the solution is 



g = V2R^\, (4.18) 

which in physical variables reads h = y^2R\z\ and thus matches to the stationary drop. 
The readily measurable amplitude ho of the thread thinning rate remains the only free 
parameter of the theory. Remarkably, the amplitude of the stress is exactly twice the 
value expc;c;ted from a naive balance of capillary and elastic forces. 

Integrating (4.15) with initial radius 1/2^/'^ given by the form of the solution (4.15) 
(plus an arbitrarily small perturbation) one finds the universal profile shape g. The 
similarity profile h is then given by (4.10) with ag 2/ij]. From Fi gure 9 one reads off 
that ho = 0.1275 for our simulation with De = 3164 ^ 1, and the numerically computed 
profiles in the corner region can easily be converted to give the rescaled profile ^(C)- 
The similarity profiles, computed from the numerical results of Figure 9 for four different 
times (solid lines), are compared with the theoretical prediction ( dotted line) in Figure 
11. As t grows, the rescaled profiles converge nicely to the theoretical prediction. 

Finally, we compare our theoretical results to experiments. To that end. high resolution 
images of experimental profiles were taken by focusing a video microscope on the corner 
region. 

Figure 12 shows a sequence of digitized images of the upper corner, which becomes 
increasingly sharp as the thread thins. To test our similarity theory, the profiles corre- 
sponding to the latest stages of pinching were rescaled using the minimum thread radius 
for both axes. As shown in Figure 13, the experimental profiles converge nicely onto 
a master curve, in very much the same way as the computed profiles of Figure 11 do. 
However, the experimental profiles turn out to be sharper ihsxi theory predicts. Thus the 
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Figure 11. A comparison of our similarity theory for largo De (dashed line) with rescaled profiles 
from a numerical simulation with parameters as in Figure 9, at logio{hmin) = —1.5, —2, —2.5, 
and —3. The dashed line was obtained from integrating (4.15), and rescaling according to (4.10), 

with = z/a^^^. 

asymptotic result given by theory needs to be shrunk horizontally by a factor of 0.5 to 
fit the experimental curves. The resulting curve is given as the dotted line in Figure 13. 

The constant ho is the only remaining adjustable parameter for the description of the 
corner region. In the case of the free jet, cf. section 2.2, it can in fact also be estimated 
in the limit that the relaxation time De of the polymer is much larger than the time tfu 
needed to form the primary filament. Namely, equating the thread thickness (3.18) given 
by the local theory of section 3 with hmin = ho exp{—tfii/3De) one finds 

ho = exp{tfu/{3De))(-j . 
Thus in the limit that is much larger than tfu one simply has 



hmin{t) = ( ) exp(-f/3£)e), 



(4.19) 
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Figure 12. An experimental closeup of the upper corner of the hquid bridge. Shown is the same 
experiment as Figure 1, but under a microscope. The time interval between two consecutive 
profiles is 0.5 sec. 

which is smaller by a factor of 2~^/^ than the result given in Bazilevsky et al. (1997) 
and Entov & Hinch (1997), based on a direct balance of elastic and surface tension forces. 
Fitting a straight line to the exponential thinning regime of Figure 5, we obtain ho = 
0.247, in excellent agreement with the theoretical prediction of (4.19), ho = (G/2)^/^ = 
0.232. 

5. Discussion 

In the present paper, confining ourselves to the simplest possible model for a dilute 
polymer solution, we have reached a rather complete (lesc;ription of the formation and 
subsequent development of the beads-on-string structure. Some details, of course, remain 
to be elaborated. 

At low viscosities or high surface tension. Oh <C 1 <C De our preliminary studies (e.g. 
Li & Fontelos (2003)) indicate that a complex succession of beads may be generated 
owing to inertial effects. Namely, fluid that rebounds from one of the large drops leads 
to capillary waves that travel on the thread, and can form one or more additional drops. 
After a sequence of drops has formed, these drops are subject to possible merging, or 
draining processes. 

A similar phenomenon of secondary drop formation has been described in Chang at 
al. (1999) and is referred to as 'recoil'. This recoil leads to the formation of a secondary 
filament between the secondary and main drops which, as claimed in Chang at al. (1999), 
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Figure 13. The last 9 profiles of Figure 12, with both cixes rescaled by the minimum thread 
radius. The dotted line corresponds to the asymptotic theory also shown in Figure 11, but with 
the horizontal ax.is shrunk by a factor of 0.5. 

exhibits a series of subsequent instabilities. However, for the parameters given in Chang 
at al. (1999) we found no indication for the formation of a secondary filament, but only 
for significantly lower viscosity or higher surface tension. Even if we do find a secondary 
filament, its presence is due to inertia, while Chang at al. (1999) claim their instability to 
be independent of inertial effects. Never did we find any indication of another instability 
leading to a third generation filament. These observations have been checked using both 
our implicit and explicit codes, making sure that the resolution was sufficient in the thin 
transition region between filament and drop. Insufficient resolution could occasionally 
lead to instabilities that vanished as resolution was increased. 

Still within the framework of the Oldroyd-B model, but at more moderate values of 
De ~ 0(1), viscous flow effects would have to be taken into account (cf. Fig. 3). As 
indicated above, the velocity in the corner region is proportional to 1/De, thus in the 
limit of large De, viscous effects drop out and only capillary and elastic forces remain. 
In the presence of viscous forces, equation (4.9) can be applied for any finite value of 
the similarity variable z. However, we have not yet succeeded in matching the solution 
described by (4.9) with the stationary drop in the general case. The reason is that whilst 
in the similarity region close to the corner (between the thread and the drop) everything 
can be described by the thread profile h alone, this is no longer possible in the transition 
region towards the stationary drop. Instead, the coupled system involving the drop shape, 
the velocity field, and the stress has to be treated. Indeed, even in the limit of large De 
where matching can be achieved using h alone, analytical computation of v and cTzz is 
far from trivial (Fontelos (2002)). 

A result of the present similarity description is that for very long times the drop 
develops into a perfect sphere, while the thread radius shrinks to zero. This means the 
profile at the point where the thread meets the drop develops into a right angle. The 
same conclusion, but using a completely different approach without making use of the 



26 



Ch. Clasen et al. 



self-similar structure of the equations, was reached in Renardy (1994). The series of 
experimental profiles shown in Fig. 12 supports this conclusion. 

Unfortunately, the fact that the slope of the profile becomes increasingly large in the 
corner means that the lubrication equations, used throughout this paper, are no longer 
valid. While it does not seem likely that this fact changes the scaling found in the present 
paper, it is clear that the precise shape of the profile will be somewhat different for the 
full equations. Full numerical simulations will also be quite challenging to perform with 
the accuracy required to resolve the corner region in the presence of strong elastic effects. 
However, one coiild use the self-similar scaling proposed in the present paper to reduce the 
dynamical 3-D equations to a time-independent equation, analogous to (4.8). In rescaled 
coordinates, the scale of the corner is of order one, and a solution of the problem could be 
attempted numerically. For the near future, the most promising prospect appears to be 
to compare with experiments with very long, well characterised polymer solutions. Note 
that even within a full 3-D theory, the correct minimum radius cannot be obtained from a 
direct balance of elastic and surface tension forces, as it was done in earlier works (Entov 
& Hinch (1997), Bazilcvsky et al. (1997)). The reason is that both contributions have 
the same sign (compare (4.3)), and thus cannot possibly balance. Instead, a matching to 
the drop has to be performed in order to determine the total tension in the thread. 

Finally, there are myriad effects associated with departures from the Oldroyd-B model, 
some of which have already been incorporated into the description of thread thinning. For 
example, even a perfectly mono-dispersc polymer solution is expected to be described 
by a spectrum of time scales, inherent in the complex relaxation dynamics of even a 
single polymer strand (deGennes (1979)). This is not expected to be a major limitation 
for the present theory, however. Bead-spring models of dilute polymer solutions indicate 
that there is a simple recursion formula relating the longest relaxation time A and higher 
order relaxation processes, of the form = A/fc™, with m = 1.5, 2 for the Zimm or Rouse 
models respectively (Bird et al. (1987)). The local rate of stretching in the cylindrical 
thread arising from the elasto capillary balance is sufficiently weak (in dimensional form 
i = 2/ (3A)), that all other modal contributions to the stress (of the form of eqs. 2.3, 2.4) 
will decay away and not contribute to the elastic stress since the convective terms will 
be of the form {2De)v' /k"^ < 1. A detailed analysis of the effect of such a spectrum of 
relaxation times was carried out in Entov & Hinch (1997), but using an incorrect form 
of the force balance, eq. (4.3). 

The effect of the finite extensibility of a polymer chain is expected to change the picture 
in Fig. 3, however, since it bounds the maximum elastic stress that can be exerted by the 
polymer chain. This is modelled theoretically by the presence of nonlinear terms in the 
constitutive equation for the polymeric stresses, which limit their growth. Thus while 
the exponential thinning described in this paper would formally lead to breakup only 
in infinite time, a real polymer thread does in fact break. Experimental observations of 
the departure from the exponential law are found for example in Bazilevsky et al. (1997) 
and Anna & McKinley (2001). The theory in Bazilevsky et al. (1997) differs from the 
conventional one in that it credits the degradation of polymers for the departure from 
the exponential law. 

A number of different nonlinear constitutive equations which bound the maximum 
polymeric stress have been proposed, and Renardy (Renardy (2001), Renardy (2002)) 
has considered the asymptotics of a number of different models. A key feature of these 
analyses is that the thread is predicted to break in finite time when the maximum elastic 
stress can no longer balance the capillary pressure j/hmin- For dilute polymer solutions in 
particular, the nonlinear form of the force-extension curve close to full extension is well- 
established both experimentally and theoretically (Li et al. (2000), Shaqfeh et al. (2003)). 
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Analysis of extensional flow of finitely extensible nonlinear elastic (FENE) dumbbells 
results in a maximum (dimensionless) polymeric stress of order azz,max ~ {2De)Gbv' 
where b = 3Nk is the finite extensibility parameter, which is proportional to the number 
of Kuhn steps or links (Nk) in the polymer chain. In the elastocapillary necking regime we 
obtain {2De)v' = 2/3 (cf. Figure 9); a naive balance thus suggests that we require Gb > 
2/3 for elastic effects to be able to grow sufficiently large to balance capillary stresses. A 
more rigorous balance is substantially more complex and requires a consideration of the 
initial polymeric stress in the filament and the additional strain accumulated by the fluid 
in the transient process of the polymeric stress growing to the saturation value <y zz.max • 
Other, more elaborate versions of the FENE model have also been developed (Lhuillier 
(2001), Ghosh et al. (2002)) which attempt to capture additional features of the internal 
dynamics of the rapid stretching process for long flexible chains. As the concentration of 
dissolved polymer chains is increased, entanglement effects also dramatically modify the 
extensional rheology of the viscoelastic fluid. However, recent experiments have shown 
that an exponential period of capillary-induced thinning, followed by finite time breakup 
still occurs (Bhattacharjee et al. (2003)). A significant benefit of the present experimental 
configuration and accompanying analytical description is precisely that the characteris- 
tics of the final breakup process sensitively depend on the nonlinear description of the 
test fluid that is used. Thus analysis of the capillary-thinning and breakup of polymer so- 
lutions provides for a unique testing ground to better understand some of the important 
nonlinear features of viscoelastic constitutive equations at large strains. 
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visit to Essen through grant SFB237. G.H.M. and C.C. would like to thank the Dupont- 
MIT Alliance and the NSF-MRI program for supporting the experimental portion of 
this work. Daniel Bonn contributed a careful reading of the manuscript and very useful 
discussions. 



Appendix A. Bounds on the velocity 

In this appendix we estimate the size of the velocity v{zo, t) at some reference point zq 
in the crossover region (b) (cf. Figure 10); to this end we integrate (2.1) from the origin 
of the self-similar region to zq: 

h^[zQ,t) 

Using the known solution (4.15) one can estimate that the main contribution to the 
integral 

(■zn 

{h')^dz 

comes from the region where h(z,t) is of order v^, £ = exp[—t/{3De)] as usual. We 
therefore introduce the rescaled profile 4'{z) by putting 

h{z,t) = Vl^{z), (A 2) 

for which (4.15) simplifies to 



28 Ch. Clasen et al. 

to leading order in £. Using (A 2) one can compute the time derivative of h to write the 
integral / as 

"'ire I (-*^+-OT.)<^- 

Integration by parts and a change of variables from z to ^ using (A 3) yield 

f 9/12 rzo/e 



3De ^ ' 3£)e 



2 2 2^2 fVWe Ci+2(t>'' 



Here we have used that according to (4.18) (j){zo/i) ~ y'2 j^ol in the crossover region 
(b). For simplicity, we also assume that the drop radius R = 1. 

Having estimated the integral we find from (A 1), using that h{0,t) ~ hoi in the neck 
region, 

zpol uUe 
as t — > oo. Prom this, the derivatives of v are easily estimated as 

dv , ^ £ , . 

^(^°'*)-T8^ ^^'^ 

and from (2.1) one finds 

. ,x -v{zo,t)hzizo,t) - htizo,t) /Ac^ 

Finally, we use (A 4), (A 5), and (A 6) to estimate the different terms of the second 
equation of (4.7). While for example £"^11 ~ £^ /De'^, the curvature term on the right of 

(4.7) can be estimated as h^K{h} /h^ ~ ^^/^/z''/^, which is always much larger in the 

transition region (b). Therefore the right hand side of (4.7) always dominates the left 
hand side in the transition region. 
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